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Abstract 

Probability distributions of the magnetic work are computed for the 2D 
Ising model by means of Monte Carlo simulations. The system is first prepared 
at equilibrium for three temperatures below, at and above the critical point. A 
magnetic field is then applied and grown linearly at different rates. Probability 
distributions of the work are stored and free energy differences computed using 
the Jarzynski equality. Consistency is checked and the dynamics of the system 
is analyzed. Free energies and dissipated works are reproduced with simple 
models. The critical exponent 8 is estimated in an usual manner. 

1 Introduction 

The Jarzynski equality [1, 2] is one of the few exact results in the context of out-of- 
cquilibrium statistical physics. This simple and elegant relation applies to systems 
initially prepared at thermal equilibrium and then driven out-of-cquilibrium by 
varying a control parameter h from say hi to hi- The probability distribution of 
the work W extracted during the experiment is related to the free-energy difference 
AF = F(h2) — F(h\) between the two equilibrium states at values h\ and hi of the 
control parameter by: 

Remarkably, this relation gives some information about the equilibrium state at the 
value hi of the control parameter even though this state has never been reached 
by the system. For a cyclic transformation, one recovers the Bochkov-Kuzovlev 
relation {e~^ w )c yc \. = 1 [3]. As recognized by Crooks, the Jarzynski equality can 
be derived from a more general fluctuation theorem [4] . It has passed experimental, 
see for instance [5, 6, 7] as well as numerical [8] tests. In both cases, the impor- 
tance of a sufficiently accurate sampling of the tail of the probability distribution 
p{W) has been emphasized. Since the pioneering experiment on DNA by Liphardt 
et al. [5], the Jarzynski relation is now widely used in chemistry and biophysics 
to estimate equilibrium free energy differences based on experiments or short-time 
out-of-cquilibrium numerical calculations. However, the way the interaction with 
the thermal bath is taken into account in the original Jarzynski's demonstration has 
been criticized [9]: the interaction is assumed to be weak enough to be neglected. It 



1 



means that the system cannot exchange heat with the bath. Equation (1) was thus 
claimed to fail since one expects that the system tends to relax to an equilibrium 
state by exchanging heat with the bath. For Jarzynski's response, see [10]. As far as 
we know, no experiment has provided any evidence of this failure yet. In this paper, 
we present a Monte Carlo study of the 2D Ising model with a Glauber dynamics. 
For such Markovian dynamics, the interaction with the bath is properly taken into 
account by the transition rates. The criticisms do not apply to this case and the 
demonstration given by Jarzynski in appendix of reference [2] is exact. 

In the context of spin models, up to now few studies have taken advantage of the 
Jarzynski equality to calculate the free energy. Results have recently been obtained 
for a single Ising spin [11] and in the mean-field approximation [12]. We present a 
Monte Carlo investigation of the two-dimensional Ising model. Numerical details of 
the calculation are presented in the first section. The three next sections correspond 
to calculations with an initial equilibrium state at different temperatures: in the 
paramagnetic phase (J3 — 0.2), the ferromagnetic phase (/? = 0.7) and at the critical 
point (/?- 0.4407). 

2 Numerical procedure 

We study the two-dimensional Ising model defined by the usual Hamiltonian 

H{{a}) = -jJ2a i a j -hJ2a i , a, = ±1 (2) 

(id) 

where the sum extends over nearest-neighbors on a square lattice. In the following, 
the exchange coupling is fixed to J = 1. Lattices with sizes from 32 x 32 to 
128 x 128 and periodic boundary conditions have been considered. The system is first 
prepared in an equilibrium state without magnetic field using the Swendsen-Wang 
cluster-algorithm [13]. In the paramagnetic and ferromagnetic phases, the system 
was first equilibrated using 200 Monte Carlo Steps (MCS). At the critical point, we 
used 1000 MCS to circumvent the critical slowing-down. This last value is a safe bet 
since it is more than two orders of magnitude larger than the autocorrelation time 
te — 5.87(1) at T c and for L — 128 [14]. We then let the system evolve according to 
a local dynamics, i.e. using the Metropolis algorithm [15], during n^er. Monte Carlo 
iterations while the magnetic field h is linearly increased. In the following, h denotes 
the final magnetic field reached after nj ter . iterations. The work is calculated as 

r-h f-tj 7 ft-i ter. — 1 

W = - Mdh = - / M(t)hdt ~ V M l (3) 

JO JO niter. i=Q 

i.e. as the average magnetization during these ni ter . iterations times the magnetic 
field h. Note that W is the mechanical work and not the work as usually defined 
in thermodynamics [16]. Let us emphasize that the demonstration of (1) given by 
Jarzynski at the end of his second paper on the subject [2] requires a well-defined 
protocol: first a sudden change Ah of the magnetic field is performed leading to a 
work —AhMi where Mi is the magnetization and second a Monte Carlo iteration 
is made which allows the system to relax and exchange heat Q = AE = E i+i — Ei 
with the bath. Note that when n;t C r. = 1 ; the Jarzynski relation is equivalent to the 
thermodynamic perturbation. The experiment is repeated n cxp = 100, 000 times. 
Instead of restarting the whole simulation, we used the last spin configuration ob- 
tained without magnetic field and did 10 further Swendsen-Wang MCS at T ^ T c 
and 50 at T c . This last value is almost ten times larger than the autocorrelation 
time. Autocorrelations between two successive experiments are thus smaller than 
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exp(-50/5.9) 
puted as 



2.10 4 and will be neglected in the following. Averages are com 
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(f(W)) = / f(W)p(W)dW = 



lox P- o=l 



f(W a ) 



(4) 



where {W Q } Q is the set of values obtained when repeating n cxp . times the numerical 
experiment. Errors are estimated as -\/[(/ 2 ) — (f) 2 ]/n cxp . as expected from the 
central limit theorem for uncorrelated random variables. 

2.1 Paramagnetic phase 

As the temperature is increased in the paramagnetic phase, the correlation length 
becomes smaller and eventually gets smaller than the lattice spacing. One can thus 
consider the system as a set of free spins. The dynamical evolution of the proba- 
bility distribution of each of them is governed by the Glauber master equation [17]. 
Since individual spins relax very rapidly to equilibrium, the probability distribution 
of the work p(W) is expected to be Gaussian [11]. It is indeed what is observed 
numerically for (3 — 0.2, as can be seen in Figure 1. We have checked that there is 
no observable deviation from a Gaussian law in the tails. Calculation of the coef- 
ficient of excess 7 = ((W — (W)) 4 ) / ((W — (W)) 2 ) 2 — 3 gives small values between 
3.10 -3 and 10 -2 for the different simulations performed. 
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Figure 1: Probability distribution of the work in the paramagnetic phase (/3 = 0.2) 
when applying linearly a magnetic field h = 0.1 at different rates: n^er. = 2,5,10 
and 20 (from right to left). The monotonous curves correspond to Gaussian fits. 



The free-energy AF can be measured using the Jarzynski equality (1), either 
directly or by using the assumption that the work W is Gaussian distributed: 



-pw\ 



1 



-pW-(W-(W)f/2^ wdw ^ AF = {w) 



2k B T 



(5) 



This last relation was first obtained by Hermans [18]. In both methods, the free 
energy difference should not depend on the experimental protocol, in our case the 
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h 




(W) 


AF (Jarzynski) 


W diss . = (W) - AF 


AF (Gaussian approx.) 


0.1 


2 


-8.83(6) 


-4.5(1). 10 1 


3.6(l).10 i 


-4.69(3). 10 1 


0.1 


5 


-2. 155(5). 10 1 


-4.70(10). 10 1 


2.5(1).10 1 


-4.70(3). 10 1 


0.1 


10 


-3.068(4). 10 1 


-4.8(1). 10 1 


1.7(1). 10 1 


-4.67(3). 10 1 


0.1 


15 


-3.523(3). 10 1 


-4.68(1). 10 1 


1.16(2).10 1 


-4.68(3). 10 1 


0.1 


20 


-3.777(3). 10 1 


-4.678(9). 10 1 


9.0(1) 


-4.68(3). 10 1 


0.1 


30 


-4.059(3). 10 1 


-4.683(5). 10 1 


6.23(8) 


-4.69(2). 10 1 


1.0 


2 


-8.578(6). 10 2 


-1.667(5). 10 3 


3.7(1).10 3 f 


-4.6(1). 10 3 


1.0 


5 


-2.0500(5). 10 3 


-2.616(2). 10 3 


2.4(2). 10 3 j 


-4.5(2). 10 3 


1.0 


10 


-2.9068(4). 10 3 


-3.523(5). 10 3 


1.5(2).10 3 f 


-4.4(2). 10 3 


1.0 


15 


-3.3037(3). 10 3 


-3.711(5). 10 3 


1.0(2). 10 3 f 


-4.3(2). 10 3 


1.0 


20 


-3.5270(3). 10 3 


-3.864(4). 10 3 


0.8(2). 10 3 f 


-4.3(2). 10 3 


1.0 


30 


-3.7663(2). 10 3 


-4.019(3). 10 3 


0.5(2).10 3 f 


-4.3(2). 10 3 


10. 


2 


-4.5445(5). 10 4 


-5.2532(5). 10 4 




-3.0(5). 10 b 


10. 


5 


-9. 1721(3). 10 4 


-9.5794(5). 10 4 




-1.9(6). 10 5 


10. 


10 


-1. 12169(2). 10 5 


-1. 14980(5). 10 5 




-1.5(5).10 5 


10. 


15 


-1. 19695(2). 10 5 


-1. 21706(4). 10 5 




-1.5(4).10 5 


10. 


20 


-1.23634(1). 10 5 


-1. 25357(5). 10 5 




-1.4(3). 10 5 


10. 


30 


-1. 27737(1). 10 5 


-1. 29189(5). 10 5 




-1.4(3). 10 5 



Table 1: Estimates of the average work (W), the free energy difference AF, the 
dissipated work (W) — AF and the free energy difference using the Gaussian ap- 
proximation in the paramagnetic phase ((3 = 0.2) for different magnetic fields h and 
transformation rates. The estimates of the dissipated work Wdiss. marked with f 
have been computed using AF as given by the Gaussian approximation because the 
Jarzynski equality fails in this case to give stable values. In the case h = 10, we do 
not even give any estimate since the Gaussian approximation leads to errors bars 
larger than the estimate of Wdiss.- 

rate at which the magnetic field is grown. Our numerical results are summarized in 
Table 1. For a small magnetic field h — 0.1 (to be compared with J = 1), the two 
methods give, as expected, estimates of AF in agreement independently of niter., 
i.e. of the rate at which the transformation is performed. Using the Jarzynski 
relation (1), the slower the transformation, i.e. the more reversible, and the more 
accurate the estimates are. Due to the amplification with the factor exp(— j3W), the 
relevant information is in the negative tail of the distribution. Insufficient sampling 
of this tail may lead to deviations of the estimate of AF. However, these deviations 
remain relatively small in our system, as can be seen in Table 2. The worst case 
is AF = —40(2) obtained for r^ter. = 2 and n oxp . = 1562 while AF ~ —46.8 for a 
larger statistics. The Gaussian assumption (5) turns out to lead to estimates of AF 
less noisy than equation (1) and that do not seem to depend on ni ter . (see Table 1). 
The relevant information comes in this case from the more probable part of the 
distribution. However, as can be seen in Table 2, the estimate gets noisier faster 
for low statistics than using the Jarzynski relation (1). We also checked that the 
free energy differences are extensive (Table 3). When intermediate magnetic fields 
h = J are applied, the Jarzynski equality fails to give estimates of AF independent 
of the transformation rate while it is still the case for the Gaussian approximation. 
When larger magnetic fields are applied, the two methods fail to give estimates of 
AF for large transformation rates, i.e. ni ter . small. A much larger statistics would 
be then required. As well-known in the context of thermodynamic perturbation, 
one can circumvent the problem by dividing the interval [0; h] in smaller pieces and 
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resorting to several simulations to estimate AF in each of them. 





(W) 


AF (Jarzynski) 


W diss . = (W) - AF 


AF (Gaussian approx.) 


1562 


-3.80(2).10 i 


-4.67(5). 10 1 


8.8(7) 


-4.7(2). 10 1 


3125 


-3.76(2). 10 1 


-4.63(3). 10 1 


8.7(5) 


-4.6(1). 10 1 


6250 


-sjstij.io 1 


-4.62(2). 10 1 


8.4(3) 


-4.6(1). 10 1 


12,500 


-3. 765(8). 10 1 


-4.63(2). 10 1 


8.7(3) 


-4.65(7). 10 1 


25,000 


-3. 773(6). 10 1 


-4.67(2). 10 1 


9.0(2) 


-4.67(5). 10 1 


50,000 


-3. 777(4). 10 1 


-4.68(1). 10 1 


9.0(2) 


-4.68(4). 10 1 


100,000 


-3.777(3). 10 1 


-4.678(9). 10 1 


9.0(1) 


-4.68(3). 10 1 



Table 2: Estimates of the average work (W), the free energy difference AF, the 
dissipated work (W) — AF and the free energy difference using the Gaussian ap- 
proximation in the paramagnetic phase (3 = 0.2 with h = 0.1 and n^er. = 20 versus 
the number of experiments n cxp . used for the calculation of the averages. 



L 


T^itcr. 


AF/L 2 (Jarzynski) 


AF/L 2 (Gaussian approx.) 


32 


2 


-2.86(2). 10~ 3 


-2.87(3). 10" 3 


64 


2 


-2.87(4). 10" 3 


-2.86(2). 10- 3 


128 


2 


-2.73(6). 10- 3 


-2.86(2). 10- 3 


32 


5 


-2.85(1). 10- 3 


-2.84(2). 10" 3 


64 


5 


-2.85(1). 10- 3 


-2.84(2). 10" 3 


128 


5 


-2.87(6). 10- 3 


-2.87(2). 10" 3 


32 


10 


-2.85(1). 10- 3 


-2.86(2). 10- 3 


64 


10 


-2.849(7). 10" 3 


-2.85(1). 10" 3 


128 


10 


-2.94(7). 10" 3 


-2.85(2). 10- 3 


32 


20 


-2.865(8). 10- 3 


-2.87(1). 10- 3 


64 


20 


-2.850(5). 10- 3 


-2.85(1). 10- 3 


128 


20 


-2.855(5). 10- 3 


-2.85(2). 10- 3 



Table 3: Estimates of the free energy difference AF in the paramagnetic phase 
((3 = 0.2 and h = 0.1) for different lattice sizes L and transformation rates h/rnt eT .- 



The results presented above can be understood by assuming that the system 
behaves as N independent Ising domains with an average magnetic moment m. 
The free energy difference between the equilibrium states with and without magnetic 
field h reads 

AF = F(h)- F(0) = -Nk B T In cosh (6) 

kbT 

The numerical data are in very good agreement with this law as can be seen in 
Figure 2. The best fit is obtained for the parameters TV ~ 0.13L 2 and m ~ 4.62 
which means that the characteristic size of the domains is of the order of £ <~ 
\J L 2 /N ~ 2.8 and that the domains arc not at saturation, Nm/L 2 ~ 0.6. At 
equilibrium, the magnetization is expected to be 
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Figure 2: On the left, free energy density AF/L 2 versus the final applied magnetic 
field h in the paramagnetic phase (J3 = 0.2) for the slowest transformation rate con- 
sidered (niter. = 30). The curve corresponds to a non-linear fit with the expression 
(6). On the right, dissipated work W^iss. with respect to the number of Monte Carlo 
iterations niter, for h = 0.1. The curve is a two-parameter non-linear fit with the 
expression (11). 

Let us assume the phenomenological evolution equation 

gM M(t)-M c , W ^ M{t) = M{0)e -t/r + l f MeM t'))e-^dt' 
at t t Jq 

(8) 

where r is the relaxation time. The work extracted from the system when the 
magnetic field is increased linearly, i.e. h(t) — ht 1 reads 



W = -h M(t')dt 



-ftM(0)(l 



-t/r\ 



h 



t r t' 



M cq .{h(t"))e- < - t '- t "V T dt'dt" 



T Jo Jo 

-hMmi-e-^)-^^- f f fe-^'-^dt'dt"^) 
rk B T J J 



h 2 Nm 2 
k B T 



t 2 -Tt + r 2 (l-e-t/r 



where M cq , has been replaced by its weak-field expansion (7). Starting the exper- 
iment from the paramagnetic phase, i.e. M(0) — 0, and subtracting the reversible 
work 

h 2 Nm 2 



Wrcv. =-h I M cq .(t')*'-- ; 



2k B T 



(10) 



one gets finally 



W dh 



h 2 Nm 2 



rt - t 2 1 - e 



-t/r 



hM cq .(h)T 



(11) 



As can be seen in Figure 2, the expression fits well the numerical data for h = 0.1. 
The non-linear fit gives r ~ 2.23 for the relaxation time and m ~ 4.72 assuming 
N = 0.13X 2 . Note that the estimate m is compatible with the one calculated from 
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AF (to ~ 4.62). The estimate of the relaxation time r is small as expected for such 
a high temperature. For larger magnetic fields, the dissipated work is still nicely 
fitted by the expression (If) but the fit gives too small average moments to. This 
is due to the fact that the expression (f f) was derived in the weak magnetic field 
limit. 

2.2 Ferromagnetic phase 

In the ferromagnetic phase (3 = 0.7, the probability distributions p(W) display two 
well-separated peaks corresponding to the work performed by systems in each one 
of the two initial ferromagnetic ground states. The average {e~@ w ) is dominated 
by the peak at more negative values of the work. As a consequence, the Gaussian 
approximation applied only to this peak leads to estimates of AF in very good 
agreement with those obtained with the Jarzynski equality (Table 4). Since the 
two peaks are not centered around zero, averages over the most negative peak were 
computed in practise by selecting configurations for which the work W is smaller 
than the average of the smallest and largest works. A large gap separating the 
two peaks, the average work would have given the same results. As can be seen in 
Table 4, the estimates of AF do not show any dependence on the transformation 
rate and give compatible values using the Jarzynski equation (1) and the Gaussian 
approximation (5). 



h 


^iter. 


(W) 


AF (Jarzynski) 


W diss . = (W) - AF 


AF (Gaussian approx.) 


0.1 


10 


1.0(5).10 i 


-1.622740(7). 10 3 


1.633(5). 10 3 


-1.624(5). 10 3 


0.1 


50 


7(5) 


-1.622732(5). 10 3 


1.629(5). 10 3 


-1.624(2). 10 3 


0.1 


250 


4(5) 


-1.622735(5). 10 3 


1.627(5). 10 3 


-1.624(1). 10 3 


1.0 


10 


-2.4(5). 10 2 


-1.63057(3). 10 4 


1.607(5). 10 4 


-1.63(4). 10 4 


f.O 


50 


-1.40(5). 10 3 


-1.630612(4). 10 4 


1.490(5). 10 4 


-1.63(2). 10 4 


f.O 


250 


-6.22(3). 10 3 


-1.6306156(9). 10 4 


1.009(3). 10 4 


-1.631(7). 10 4 


to. 


fO 


-1.076(2). 10 b 


-1.63676(1). 10 b 


5.61(2). 10 4 


-1.6(2). 10 b 


to. 


50 


-1.3313(10). 10 5 


-1.637355(8). 10 5 


3.061(10). 10 4 


-1.64(7). 10 5 


to. 


250 


-1.4475(6). 10 5 


-1.637402(2). 10 5 


1.899(6). 10 4 


-1.64(3). 10 5 



Table 4: Estimates of the average work (W), the free energy difference AF, the 
dissipated work (W) — AF and the free energy difference using the Gaussian ap- 
proximation with the most negative peak in the ferromagnetic phase (/3 = 0.7) for 
different magnetic fields h and transformation rates h/ niter.- 



In a sense, the physics is simpler than in the paramagnetic phase. Before apply- 
ing the magnetic field, the initial magnetization is either +M or — M with small 
fluctuations around these two values. When applying a small magnetic field, i.e. 
h = 0.1 in our case, the magnetization does not change much so that the work is 
either —M h or +M h. According to Jarzynski equation (1), the free energy reads 



AF = -k B T\n{e- fiw ) = -k B T\n 



}_ e -f3M a h + l e /3M h 



~ -M h+fc B Tln2 (12) 



Using the saturation magnetization M sat = L 2 = 16384, one obtains indeed a good 
estimate of AF. The average work is small (W) = (—M n h + M h)/2 = so that the 
dissipated work is Wdiss. = (W) — AF ~ — AF. For larger magnetic fields, h — 1 and 
h = 10, this approximation is not valid anymore: the magnetization tends to reverse 
itself if it was initially in the opposite direction of the field. However, since the 
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Jarzynski equality (1) is dominated by the systems for which the magnetization was 
initially in the same direction than the field, one still expects that AF ~ — M sat /i. 
But the average magnetization now moves towards — M sat so Wdi ss . > 0. Obviously 
dissipation is thus mainly due to the systems for which the initial magnetization 
was in the "wrong" direction. 

2.3 Critical point 

We now present results at the critical point (3 C = |ln(l + \/2) — 0.440687. Like 
in the ferromagnetic phase, the probability distributions p(14 7 ) display two well- 
separated peaks. Due to finite-size effects, the initial state has indeed a non- 
vanishing magnetization of sign either positive or negative due to the Z2 symmetry 
(<Tj — ► — <jj) of the Ising model. In contradistinction to the ferromagnetic case, the 
two peaks do not have a Gaussian shape. The free energy differences AF calculated 
from the Jarzynski equation (1) are presented for two different magnetic fields in 
Table 5. For the smallest magnetic field h ~ 0.023, the free-energy turns out to be 
independent of the transformation rate h/riitcv.' the maximum deviation from the 
average is of the order of one standard deviation. For the largest magnetic field 
h ~ 0.23, a systematic deviation is observed. 



h 


^iter. 


(W) 


AF (Jarzynski) 


W diss . = (W) - AF 


0.022692 


2 


-1.6(7) 


-2.690(7) 


TO 2 


2.67(1).10 2 


0.022692 


10 


-1.1(7) 


-2.697(4) 


.10 2 


2.69(1).10 2 


0.022692 


50 


-7.5(7) 


-2.704(5) 


.10 2 


2.63(1).10 2 


0.022692 


100 


-1.39(7). 10 1 


-2.701(4) 


.10 2 


2.56(1).10 2 


0.022692 


150 


-1.85(7). 10 1 


-2.698(2) 


.10 2 


2. 513(8). 10 2 


0.022692 


250 


-2.87(7). 10 1 


-2.699(1) 


TO 2 


2.411(8). 10 2 


0.226918 


2 


-2.4(7). 10 1 


-2.897(2) 


.10 3 


2. 873(9). 10 3 


0.226918 


10 


-1.87(7). 10 2 


-3.009(2) 


.10 3 


2. 822(9). 10 3 


0.226918 


50 


-7.93(6). 10 2 


-3.103(1) 


.10 3 


2. 310(7). 10 3 


0.226918 


100 


-1.362(5). 10 3 


-3.133(2) 


.10 3 


1.771(6). 10 3 


0.226918 


150 


-1.680(4). 10 3 


-3.136(2) 


.10 3 


1.456(6). 10 3 


0.226918 


250 


-2.017(3). 10 3 


-3.143(1) 


.10 3 


1.126(4). 10 3 



Table 5: Estimates of the average work (W), the free energy difference AF, the 
dissipated work (W) — AF at the critical point for two different transformation 
rates h/mt eT . and two magnetic fields h. 

We made simulations for ten values of the final magnetic field in the range 
h € [0.023; 0.23]. The singular part of the free-energy is expected to scale at the 
critical point as 

F sing .(h) ~ h 1+1 ' s (13) 

T->T C , 

Neglecting the contribution of the regular part, we fitted AF as given by the Jarzyn- 
ski equality with the scaling form (13). The data display a nice power-law behavior 
but fluctuations are observed for fast transformations (niter. < 10) at large magnetic 
fields. The estimates of 1 + 1/5 and S are collected in Table 6. As the transformation 
gets slower, and thus niter, larger, the numerical estimate gets closer to the exact 
result (5 — 15. For the fastest transformation ni ter . = 2, the relative deviation of the 
numerical estimate from the exact result is of order of 100%! Note however that 
much better estimates are obtained when restricting the fit to the smallest values 
of the magnetic field: 6 = 19.9(7) for h < 0.018 and S = 14.1(6) for h < 0.013 when 
niter. = 2. Since error bars take into account all sources of statistical fluctuations, 
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the deviation may be explained as a systematic bias due to a too small number of 
experiments for the calculation of the average. As can be seen in Table 7, when n oxp . 
is made smaller the estimate of 6 indeed increases. The effect is dramatic for fast 
transformations: a relatively good estimate of 6 is already obtained for n^er. = 250 
when averaging over only n oxp . = 1000 experiments but the same deviation from 
the exact result 6 = 15 is not even obtained with n oxp . = 100,000 experiments for 
niter. — 50. More stable estimates of 6 are obtained for a smaller system, L = 64, 
for which 6 does not display a systematic deviation but fluctuates between 14.0(1) 
{nuer = 2) and 14.1(2) (n iter = 150) and remains below the exact value because of 
finite-size effects. 





1 + 1/5 


6 


2 


1.024(7) 


42(1) 


10 


1.045(6) 


22.2(3) 


50 


1.059(5) 


17.0(1) 


100 


1.064(4) 


15.7(1) 


150 


1.065(3) 


15.36(7) 


250 


1.066(2) 


15.13(5) 



Table 6: Critical exponents 1 + 1/6 and 6 obtained by power-law interpolation of 
the free energy AF = F(h) — F(0) <~ h 1+1/>s with ten magnetic fields in the range 
[0.023; 0.23] at the critical point for different transformation rates h/nnci.- 





1 + 1/5 


6 


10 2 


1.060(10) 


16.6(3) 


10 3 


1.066(7) 


15.1(2) 


10 4 


1.065(3) 


15.44(8) 


10 5 


1.066(2) 


15.13(5) 



Table 7: Critical exponents 1 + 1/6 and 6 obtained by power-law interpolation of 
the free energy AF = F(h) — F(0) ~ h 1+1 / s with ten magnetic fields in the range 
[0.023; 0.23] at the critical point for different number of experiments n oxp . in the 

slowest ('MSG ^Ijtcr 

. = 250. 

As already mentioned, the two-peak structure of the probability distribution of 
the work p(VF) differs from the ferromagnetic case. As can be seen in Figure 3, the 
left peak, which dominates the average in the Jarzynski equality (1) and thus the 
free energy difference AF, moves monotonously to the left as the magnetic field is 
increased. The right peak moves first to the right, up to a certain magnetic field 
before pursing to the left. This behavior can be understood by assuming that the 
dynamics of the magnetization is governed by the Langevin equation 

dM 6T /■* 

^t = -jm^ Kh(t) * M[t) - M(0) + K I mdt ' (14) 

where the link with equation (8) is made by considering k as a magnetic suscepti- 
bility divided by the relaxation time r. It is sufficient for the discussion to consider 
k as constant. For a single realization, the initial magnetization M(0) is non-zero 
due to finite-size effects. When the field is increased linearly, i.e. h(t) — fit, the 
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work extracted from the system up to time t is 

W = - [ M(t)hdt = -M(0)ht - ^h 2 t 3 = -M(0)h - ^h 2 t (15) 
Jo 6 6 

When M(0) < 0, the work initially increases and eventually at h t = 3\M(0)\/nt 
starts to decrease as can be observed in Figures 3 and 4. In our case, no turning point 
is observed for niter. = 10 because the magnetic fields remain too small. The data 
presented in the inset of figure 4 reproduce the linear behavior of (W/h — Wo)/ni te r. 
versus h predicted by equation (15). Higher order terms are observed only for the 
slowest transformations. When the initial magnetization is initially in the same 
direction as the magnetic field, i.e. M(0) > 0, the work is always negative. For 
very slow transformation rate, the magnetization is expected to follow the magnetic 
field, i.e. to behave as the equilibrium magnetization M cq . <~ h x l & . Taking into 
account the magnetic field dependence of the k, one should expect in this case the 
solution of the Langevin equation (14) to scale as 

Wrcv. ~ - / h' 1/S dh' ~ -h 1+ * (16) 
Jo 

This behavior is indeed observed: a power-law interpolation of the position of the 
left peak versus the magnetic field leads to exponents S = 14(1) for rji t0 r = 150 and 
250 compatible with the exact result 6 = 15. 




Figure 3: On the left, probability distributions p(W) versus the work W for 
a fast transformation niter. = 10 for four different final magnetic fields h — 
0.004,0.017,0.031 and 0.044 (from bottom to top). On the right, probability dis- 
tribution p(W) versus W for an intermediate transformation ni te r. = 100 for the 
same values of the magnetic field (from bottom to top) . 



3 Conclusions 

The Jarzynski equality offers the possibility to estimate free energy differences in a 
very simple manner in Monte Carlo simulations. Despite the fact that this relation 
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I 1 1 1 1 1 1 1 1 1 1 ^000 

0.05 0.1 0.15 0.2 0.25 

h 

Figure 4: Average work (W) of the right peak of the probability distribution at the 
critical point versus the final magnetic field h for different transformation rates, i.e. 
niter, from 2 to 250. In the inset, ((W)/h — M )/n iter . where M has been adjusted 
to the value M — 9016 is plotted versus h to check the linear behavior predicted 
by (15). 

is exact in the case of a Markovian dynamics, it has to be used carefully. As already 
shown by several authors, we observed systematic deviations for strongly irreversible 
transformations due to an insufficient sampling. These systematic deviations are 
caused by the fact that the average (e~@ w ) is dominated by very rare events — W 3> 
1 that may require a large number of simulations to be properly sampled. This 
limitation of the Jarzynski equality is a strong one: the error cannot be estimated 
simply. When the distribution of the work is Gaussian, it has been shown that 
the number of experiments has to be grown exponentially with the dissipated work 
Wdiss. [19]. We have observed that the Gaussian approximation (5) may still give 
reliable estimates of the free energy difference AF when the Jarzynski equality fails 
as for example in the case of intermediate magnetic fields applied on a paramagnet. 
Unfortunately, the Gaussian approximation is valid only for one-particle systems [20, 
21] or systems with short-range correlations like paramagnets but not at the critical 
point. However, the Jarzynski equality turned out to be useful at small magnetic 
field. It was reliable enough for instance to estimate the critical exponent 8. For this 
purpose, the Jarzynski equality may be useful since a single Monte Carlo simulation 
is required if the work is measured at different times. 
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